﻿using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Linq;
using System.Numerics;
using System.Runtime.Intrinsics;
using System.Runtime.Intrinsics.X86;
using System.Text;
using System.Threading.Tasks;

namespace LinearFunctionFitting.Demo.Algorithm
{
    public unsafe class VectorMethods : IFunction
    {
        public (int a, int b) LinearFunction(int[] xs, int[] ys)
        {
            int n = xs.Length;
            int a = 0;
            int b = 0;
            int vsize = Vector256<int>.Count;
            fixed (int* xptr = xs, yptr = ys)
            {
                Vector256<int> vSumX = new Vector256<int>();
                Vector256<int> vSumY = new Vector256<int>();
                Vector256<int> vSumXY = new Vector256<int>();
                Vector256<int> vSumX2 = new Vector256<int>();
                int i;
                int limit = n - (n & (vsize - 1));
                for (i = 0; i < limit; i += vsize)
                {
                    Vector256<int> vx = *(Vector256<int>*)(xptr + i);
                    Vector256<int> vy = *(Vector256<int>*)(yptr + i);
                    vSumX = Avx2.Add(vSumX, vx);
                    vSumY = Avx2.Add(vSumY, vy);
                    vSumXY = Avx2.Add(vSumXY, vx * vy);
                    vSumX2 = Avx2.Add(vSumX2, vx * vx);
                }
                int sumX = 0;
                int sumY = 0;
                int sumXY = 0;
                int sumX2 = 0;
                int x2Sum = 0;

                for (; i < n; i++)
                {
                    int xt = *(xptr + i);
                    int yt = *(yptr + i);
                    sumX += xt;
                    sumY += yt;
                    sumXY += yt*xt;
                    sumX2 += xt*xt;
                }
                int* psumX = (int*)&vSumX;
                int* psumY = (int*)&vSumY;
                int* psumXY = (int*)&vSumXY;
                int* psumX2 = (int*)&vSumX2;
                for ( i = 0; i < vsize; i++)
                {
                    sumX += *(psumX + i);
                    sumY += *(psumY + i);
                    sumXY += *(psumXY + i);
                    sumX2 += *(psumX2 + i);
                }
                x2Sum = sumX * sumX;
                int nsumX2 = n * sumX2;
                int k = x2Sum - nsumX2;
                b = (sumX * sumXY - sumY * sumX2) / k;
                a = (sumY * sumX - n * sumXY) / k;
                return (a, b);
            }
        }
    }
}
